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Abstract. An instability driven by the thermal anisotropy of a single electron species 
is investigated in a 2D particlc-in-cell (PIC) simulation. This instability is the one 
considered by Weibel and it differs from the beam driven filamentation instability. A 
comparison of the simulation results with analytic theory provides similar exponential 
growth rates of the magnetic field during the linear growth phase of the instability. 
We observe in accordance with previous works the growth of electric fields during 
the saturation phase of the instability. Some components of this electric field are not 
accounted for by the linearized theory. A single-fluid-based theory is used to determine 
the source of this nonlinear electric field. It is demonstrated that the magnetic stress 
tensor, which vanishes in a ID geometry, is more important in this 2-dimensional 
model used here. The electric field grows to an amplitude, which yields a force on 
the electrons that is comparable to the magnetic one. The peak energy density of 
each magnetic field component in the simulation plane agrees with previous estimates. 
Eddy currents develop, which let the amplitude of the third magnetic field component 
grow, which is not observed in a ID simulation. 
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1. Introduction 

Plasma instabilities that result in the growth of magnetic fields from noise or in their 
amplification are important for the magnetization of the upstream region of supernova 
remnant shocks [1] or for the cosmological magnetic field generation [2]. Two kinetic 
instabilities, which are based on idealized electron distributions, are frequently discussed 
in the context of providing the seed fields for further instabilities. These are the Weibel 
instability and the filamentation (beam- Weibel) instability. 

In a plasma with a thermally anisotropic electron distribution, in which the electron 
thermal spread in one direction is larger than that in the other two, electromagnetic 
fluctuations (noise) are amplified by the thermal anisotropy-driven Weibel instability 
(TAWI). A magnetic field is generated, which is coherent on electron skin depth scales. 
This plasma system is stable against electrostatic instabilities, as long as there is no 
drift between the electrons and ions. The filamentation (FI) or beam- Weibel instability 
[3-9] acts in systems with counter-streaming electron beams and it can be considered as 
an extreme form of the TAWI. The FI is efficient in plasmas with relativistic streaming 
velocities between electron beams with a comparable density, because the growth time of 
the electromagnetic instability is small compared to that of the electrostatic instabilities. 
Recently, the importance of the TAWI and the FI has also been recognized in the fast 
ignition processes for the inertial confinement fusion [10]. 

Instabilities driven by a thermal anisotropy have been widely examined numerically 
and analytically [11-17] and they are considered also in our present work. The 
thermal velocity of the bi-Maxwellian distribution is larger along one axis than in its 
perpendicular plane. This anisotropy induces higher micro-currents along the hotter 
direction and their magnetic repulsion separates in space the electrons with oppositely 
directed velocity vectors. The low thermal energy in the perpendicular plane cannot 
work against the structure formation in this plane. A net current and electromagnetic 
fields develop. The magnetic energy density can reach in extreme cases up to 1/12 of the 
total energy density [12-14] and it exceeds by far the electric one. The focus of previous 
studies has thus been on the magnetic field. The neglect of the electric field is probably 
justified, if the thermodynamic properties of the plasma are considered. However, if 
the electron speeds are well below the speed of light, the electric forces on individual 
electrons may not be small compared to the magnetic ones. This has motivated several 
recent investigations of the nonlinearly driven electric field. Vlasov and PIC simulations 
[15-17] have demonstrated, that the electric and magnetic field structures are linked. 
It turns out that the driver of the electric field is the pressure gradient force of the 
self-generated magnetic field, if the wave spectrum is limited to one dimension [17]. 

The magnetic pressure gradient force and the magnetic tension force are both a 
consequence of the J x B force of the self-generated magnetic B-field on the driving 
current J. Only the magnetic pressure gradient force can, however, develop in the system 
with a one-dimensional wave spectrum investigated in Ref. [17]. The magnetic pressure 
gradient force remains stronger than that due to the magnetic tension as we go from 
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ID to 2D simulations of the FI [18], but it is unknown if this is also true for the TAWI. 
We address this issue here. We consider immobile ions and a bi-Maxwellian electron 
distribution with a large temperature along one axis (the parallel component in the 
following) and a lower temperature in the perpendicular plane. The plasma parameters 
are similar to those in Ref. [17], but here the simulation geometry gives rise to a two- 
dimensional wave spectrum and, thus, to a magnetic tension force. We can determine 
the relative importance of both components of the J x B-force for the electric field 
generation by their direct comparison. We may also expect consequences of the altered 
filament dynamics in a 2D simulation. The magnetic field in a ID geometry eventually 
becomes strong enough to keep filaments separated, suppressing their further merging. 
A 2D geometry allows repelling filaments (oppositely directed current) to move around 
each other and continue to merge with attractive filaments. 

Often, the magnetic trapping mechanism is invoked to explain the saturation of the 
instability, i. e. the condition for saturation is given, when the magnetic bounce frequency 
wb = IqkvBk/mcl 1 / 2 is of the same order as the growth rate of the TAWI [23]. This 
magnetic trapping mechanism does, however, not take into account the electric field. It 
is, however, becoming increasingly evident that the electric forces can not be neglected, 
when the TAWI or the FI saturate [24]. The pressure gradient of the magnetic field 
driven by the FI of counter-propagating electron beams accounts for the electric field 
in ID and 2D simulations [18, 25]. This electric field is driven by the magnetic pressure 
gradient force, if the wave spectrum of the TAWI is one-dimensional [17]. Early 2D PIC 
simulations of the TAWI [11] did not have the signal-to- noise ratio that is necessary to 
determine the exact source of the electric field [11] and the recent simulations by [15- 
17] did not resolve the 2D wave spectrum of the TAWI. It is thus not well-understood, 
which mechanism produces the electric fields if the wave spectrum of the TAWI is two- 
dimensional. For this purpose the spatio-temporal evolution equation of a single-species 
fluid is considered and the coupling between the perpendicular and parallel components 
of the magnetic and electric field is analyzed. 

The linear theory and the numerical method are outlined in section 2 and the 
simulation parameters are specified. In section 3 the numerical results are presented, 
where we get a similar power spectrum to that observed for the FI [19]. The magnetic 
energy density exceeds by far the electric energy density in our simulation. However, 
the magnitudes of the electric and magnetic forces, which act on the non-relativistic 
electrons, are comparable; both are thus equally important for the particle dynamics. 
Most importantly, we find that the magnetic pressure gradient force by itself is too 
weak in the 2D simulation to explain the observed electric field. The electric force on 
the nonrelativistic plasma electrons is instead comparable to that expected from the 
superposition of the magnetic tension force and the magnetic pressure gradient force. 
The spatial distributions of the fields are alike but not identical. In section 4 the results 
and their implications are discussed. 
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2. The instability, initial conditions and the simulation method 

2.1. The linear instability 

The Weibel instability is investigated in a homogeneous, collisionless plasma with 
the initial magnetic and electric field strengths B = E = 0. The ions form an 
immobile background that compensates the electron charge. The spatially uniform 
initial distribution of the electrons is given by 



fo(v±,v\\) = (0 s 3/ 2 2 ex P 



9 2 + o 2 
zu th\\ 



(1) 



(2) 



where v t h\\ = \/kT\\/m and Vth± = \J kT± / m denote the thermal velocities of the parallel 
and the perpendicular components, respectively. The Boltzmann constant is k, m is the 
electron mass and Tj_ and Ty are the respective temperatures. 

An instability is driven by a temperature anisotropy A = (v t h\\/v t h±) 2 — 1^0. 
We choose the large A = 399 here and in the simulation (v th \\ = 20v th ±), so that 
we get a large growth rate of the instability and a good signal-to-noise ratio for the 
electromagnetic fields in the simulation. This value of A furthermore allows us to test 
for a large A the finding [12, 14], that the maximum average magnetic energy density is 
almost independent of the initial anisotropy for A > 25. Electromagnetic fluctuations 
with a wavevector in the perpendicular plane are amplified in this case, according to 
the well-known dispersion relation for the linear phase of the instability [21] 
k 2 c 2 a 2 \ 1 x ./ io \ 

uj; id 2 |_ 2 \kv th ±J 
k and a are the wave number and the associated linear growth rate of the growing 
electromagnetic oscillations with a purely imaginary frequency u = io . The electron 
plasma frequency is given by co p = {e 2 n/ e^m) 1 / 2 with electron charge e and electron 
number density n and Z'(() = — 2[1 + C-^(C)] is the first derivative of the plasma 
dispersion function Z(() = ix~ x l 2 dt exp(—t 2 )/(t — Q. The normalised wave number 
k max c/tUp = A x l 2 determines the upper limit of unstable wave numbers < k < k max . 

2.2. Numerical method and code resolution 

The PIC method models self-consistently the interplay of the electric and magnetic 

fields with a collision-less kinetic plasma. The plasma is treated as an incompressible 

phase space fluid, which is approximated by an ensemble of computational particles 

(CPs). Each CP has the same charge to mass ratio q cp /m cp as the physical particles it 

represents. With the relativistic momentum p cp = m cv ^ cv and velocity v cp of a CP, 

the Maxwell equations for the electric and magnetic fields E and B 

<9B 1 <9E 

VxE=--. VxB = -- + „ J, (3) 

and the Lorentz equation for the CP that is located at the position x cp 

= q cp (E(x cp ) + v cp x B(x cp )) (4) 
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Figure 1. The growth rate a(kc/wp)/u p as a function of the wavenumber k in the 
perpendicular plane. After a steep rise at low wavenumbers the maximum is reached. 
Thermal effects limit the growth rate at large k. The cut-off is given by A 1 / 2 « 20. 

are solved. The code fulfills V ■ E = p/e and V • B = to round-off precision [22]. 

In contrast to the freely moving CPs, the electric and magnetic fields are defined on 
a grid and have to be interpolated to the position of each CP. With Eq. (4) the velocity 
v cp is updated and the particle position is advanced in time with dx cp / dt = v cp and the 
time step A t . The total current J contains the contribution of all microcurrents g cp v C p, 
which are then interpolated back onto the grid. Then the electric and magnetic fields 
are updated with (3) and the individual steps are repeated. 

2.3. Initial conditions and simulation setup 

The parallel axis is aligned with the (unresolved) z direction {yn = v z ) of the simulation. 
According to the linear theory, the wave vectors in the x-y plane, which we label the 
perpendicular plane, are unstable and B x and B y grow exponentially. The rotational 
symmetry around the parallel direction allows us to combine them to B± = B x + iB y . 
This is allowed if the plasma structures are large compared to a grid cell and small 
compared to the box, which are both rectangular. The electric field driven through the 
displacement current is E\\ — E z . In the following, we will discuss the magnetic field in 
terms of cB± and cB\\ as they have the same unit as the electric fields E± and E\\. 

The boundary conditions are periodic in all directions. We set the electron plasma 
frequency to u p = 6.3 • 10 5 s _1 . The simulation box is composed of N g x N g = I600 2 
rectangular grid cells, each with a side length A x = 10 m or A x u p /c = 0.021 in terms 
of the electron skin depth, in the x and y directions. The total box size is therefore 
(N g A x uj p /cy = 33.5 2 . The number of CPs per cell is initially N = 160. The simulation 
time step is u p A t w 0.0094 and the simulation time is u p T sim = 434. The thermal 
velocities are v t h\\/c 9 x 10~ 3 and v t h\\ = 20v t h±- 
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3. Numerical results 

3.1. The energy densities 

An investigation of the energy densities provides the overall temporal evolution of the 
instability. Figure 2 shows the parallel and perpendicular components of the electric and 
magnetic energy densities, which are given by e£ S (£) = N~ 2 ^o\E 8 (JA x , kA x , t)\ 2 /2 

and €Bs{t) = N~ 2 ^2 \B s (jA x , kA x , t)\ 2 /2fi . The magnetic and electric fields are B x + 

iB y and E x + iE y for s =_L and B z and E z for s =\\. The energy densities are normalised 
by exo = e K {t = 0). The kinetic energy density e#(£) = N~ 2 A~ 3 Yl, m cpC 2 {.lj — 1)- 

3 




Figure 2. (Colour online) (a) The normalised magnetic energy densities €B±/f Ko and 
£b\\/cko (upper black curves) and the electric energy densities e_Ej_/ e ft"o and CE\\/ e K0 
(lower blue dashed curves). Initially only es± grows exponentially and the onset of 
its saturation at tui p « 60 triggers the growth of the other components. The energy 
density of B± peaks above a value 1/12 and just below 1/6 (horizonal line), (b) The 
decrease of the averaged anisotropy < A >. (c) The kinetic energy density of particles 
moving in parallel and perpendicular direction. 

The growth of e^j_ is exponential for tu p < 60 and we notice a weak growth of 
em, while the other energy densities are at noise levels. This defines the linear phase of 
the instability. The e^j_ ^> €e\\ at this early time, because the latter is electromagnetic 
in the chosen geometry, while €e± contains also electrostatic noise. After tu p ~ 60 the 
growth of €b± slows down, while the other energy densities grow rapidly. Processes must 
be at work, which are not captured by the solution of the linear dispersion relation. The 
6b\\ is amplified here, but not in a ID simulation that evidences only a growth of €e± 
and of €e\\ [17]. All energy densities saturate at the same time at tu p ~ 80. Furthermore, 
the energy density €e± is two to three orders of magnitude smaller than €b_l, thus the 
amplitude of the mean electric field is one to two orders of magnitude lower than the 
amplitude of the mean magnetic field. Most electrons have velocities well below c/10 
and they are therefore affected by the electric force as much as by the magnetic force. 
Neglecting the electric field [11] is thus not necessarily permitted. 
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We can approximate with an exponential function the early stage of the evolution 
of the energy densities of the perpendicular components (not shown), yielding the 
growth rates cb± ~ 0.07u p , which is below the analytical value a Ana = 0.114u p , while 
cte± ~ 0.06u p . It is the same situation as in reference [17]: The growth rate of the 
magnetic energy density is reduced in comparison to the analytical value. This is at 
least partially due to the averaging over all wave numbers. The growth rate of €e± 
is here below that of €b±, while it grew twice as fast in a ID simulation [17]. Such 
a reduction of the growth rate of e^± as we go from a ID to a 2D simulation is also 
observed for the FI [18]. 

In Figure 2 (b) the decrease of the thermal anisotropy is shown. With the transition 
to the nonlinear phase the thermal temperatures have almost equalized. Also the kinetic 
energies of particles moving in parallel and perpendicular direction balances (see Figure 
2(c)). 



3.2. The power spectra of the field components 

In what follows, we analyse the nonlinear evolution of the fields. Figures 3(a-c) display 
the power spectra of cB±, cBu and E±. The power spectra as a function of the scalar 
wavenumber k = ik\ + ky) 1 ^ 2 have been obtained by the integration over the azimuth 
angle of the power spectrum in the k x — k y plane. We have normalised the spectra by 
Pa, which is the highest value of the power spectrum P(k,t) of B±. 
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Figure 3. (Colour online) The power spectra P(/e,i)/P of cB± (a), cBu (b) and 
Ej_ (c) as a function of the logarithmic (base 10) wave number and time. Pq is the 
maximum value of the power spectrum of cB±. The colourbar is the 10- logarithm 
of the normalised power. The overplotted curve k' cx t _1 follows the evolution of the 
power maximum to low k. 



The cB± grows exponentially and to a significant power over a wide range of wave 
numbers up to a k = k maXt si m with log 10 (k maX) si m c/ujp) ~ 1.2. This is close to the 
analytical value \og 10 (k max ^Anac/to p ) = log 10 \/C4 ~ 1.3. The P/Po gradually goes over 
into noise at higher k. The magnetic field components cB± and cB» grow in a similar 
range of wave numbers, which suggests that they are coupled. This mechanism must 
be nonlinear, because the solution of the linear dispersion relation does not predict the 
growth of cBh . The correlation of both magnetic components is reduced at later times, 
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when the spectrum of cB± broadens. The wavenumber spectra of E± and of B\\ evidence 
the onset of the growth of structures at about the same time tu p ~ 10 2 . However, the 
power spectrum of the electric field reaches larger wave numbers. The electric field thus 
changes more rapidly in space than the magnetic one, as it would be the case if it were 
driven for example by the magnetic pressure gradient oc VB 2 . 

Figure 3 furthermore shows that the power of the field components during the non- 
linear phase shifts in time to lower wave numbers k. This goes along with an increase 
of the scale size (coherence length) of the magnetic field by the merging of the current 
structures in position space. The merging requires two spatial dimensions orthogonal to 
the parallel axis [11]. In a ID simulation the merging of the current filament eventually 
stalls [11, 16, 17]. The wave number associated with the maximum power of cB± can be 
approximated by the curve k' oc However, the power spectra are broadband and 
the sampling time too short to determine accurately this dependence. 




log 10 fcc/co p ■ 

Figure 4. 10-logarithmic plot of the power spectrum at tuj p — 434. The dashed lines 
are the fits to the simulation data. 

Figure 4 is a double logarithmic plot of the power spectra of the field components 
cB±, cB\\ and E± against the normalised wave number at the end of the simulation. 
Only the range up to log 10 (/cc/cu p ) = 1.15 is shown here as the noise dominates at large 
values. The power spectra have a power-law behaviour for almost the whole range. The 
indices are —5.2 for cB±, —5.4 for cB\\ and —0.9 for E±. The power spectra of B± and 
B\\ are practically identical, except for a scaling factor 10. This is further evidence 
for a nonlinear coupling between both components. In contrast, the power spectrum 
of the electric field is flat compared to the magnetic components. The amplitudes for 
cB±, cB\\ and E± are proportional to a force. The spectral distributions reveal a trend, 
namely that the electric forces can become more important than the magnetic ones, if 
the dynamics of non-relativistic electrons on small scales (large k) is considered. 
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3.3. The saturation mechanism 

The field data are now analyzed in order to give information about the saturation 
mechanism and about how the driving field B± couples to E± and to B\\. We consider 
for this purpose the equation 

1 * pry 

d t (nv) + V • (nvv) = V • P T (E + v x B), (5) 

m e m 

with which we can investigate the interplay of the field components. It describes the 
spatio-temporal evolution of a single-species fluid consisting of electrons with particle 
density n, mean velocity v and thermal pressure tensor Pp. With Ampere's law 
J = ^VxB - e ff , J = -env and (V x B) x B = - [VB 2 /2 - V • (BB)] the 
equation can be rewritten as 

d t (nv) + V-(nvv) = -— V-P T -— E--V-<x--VPij + -Bx ^(6) 

m e mm m m at 

where a = BB//i and Pb = B 2 /2/i describe the magnetic stress tensor and pressure. 

If we consider a configuration where the gradients d y and d z are not resolved, the 
term —V • a/m due to the magnetic tension vanishes. The development of the FI or 
the TAWI imply, that the term m^VPe on the right hand side of Eq. 6 becomes 
important. The PIC simulations in [17, 18] demonstrate that an electric field grows 
with an amplitude that cancels the term m~ 1 VPe by the second term enE/m on the 
right hand side of Eq. 6, provided that immobile ions are considered. The amplitude of 
the electric field along the ID simulation box is then given by 

E x = -- — B 2 = — P B 7 

2en e no dx en e dx 

As a result, the electric field energy density grows twice as fast as the magnetic one. 
The TAWI and FI saturate, when the combined electric and magnetic force prevents a 
further spatial re-arrangement of the current by confining the electrons in space. 

The perpendicular gradients are resolved (d x ^ 0, d y ^ 0) in our 2D PIC simulation 
and V • a no longer vanishes. If we neglect the contribution of the thermal pressure 
gradient and the displacement current (term 1 and 5 on the right hand side in Eq. 6), 
then the right hand side of Eq. 6 vanishes, provided that 

E ± = —V xy ■ a - —V xy P B , (8) 

en en 

where V ' xy is the gradient in the x-y plane. We will examine here the magnitude of the 
three terms in Eq. 8 using the PIC simulation data. 

Movie 1 shows the development of the individual contributions: The upper left 
panel shows the magnetic field |B(x, y)\ normalised by u p m/e. The middle panel shows 
the development of the nonlinear terms in Eq. 8 in the same normalisation oj p m/ce 
as the electric field, i. e. the superposition of the magnetic pressure gradient VPs an d 
the divergence of the magnetic stress tensor V • a. The bottom panel is the normalised 
electric field component E±(x,y). Besides the formation of larger structures in the 
magnetic field |B(x, y)\, the connection of the magnetic tension and pressure gradient 
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force with the perpendicular electric field is obvious at various locations in Movie 1. For 
tu p > 50 structures in the electric field are contrasted with the noise background and the 
strongest ones correspond to the two non-linear terms we examine here. Figure 5 displays 



-3 




X-position X-position 

Figure 5. (Colour online) The relevant field and force components at u> p t = 75 in a 
subinterval of the box: (a) \B±\ = \B X + iB y \, (b) the parallel magnetic component 
|B||| = \B Z \, (c) the \E X \ = \E x + iE y \. The force |V-<7 + VP s | (d) and the individual 
components (e) |V • a\ and (f) |VPg|. Values are normalized and can be compared 
directly. 

the most relevant field and force components at the time u p t = 75, when the TAWI has 
just saturated in Fig. 2. The dominant field component is clearly \B±\, as expected 
from the TAWI and from the energy density diagram. The peak amplitude ~ 0.02 
is compatible with a saturation by the magnetic trapping mechanism. The magnetic 
bounce frequency in the normalized units k = kc/u p , v = v/c and B± = mBj_/e 

1/2 

displayed in Fig. 5 is (ujb/u p ) = \kvB±\ . We rearrange the equation and get 
v = (ub/u p ) 2 /{\B±\k). Equating the ub with the typical growth rate from Fig. 1 
gives ub/^p = 0.05 and we further assume that k = 2.5 and \B±\ = 0.02. We get the 
reasonable rough estimate v ~ 0.05, which is a few times v t u\/c ~ 0.01. 

The correspondence between the total force |V-<7+VPb| and \E±\ is apparent. The 
amplitudes are very similar and at least the strong structures in Figs. 5(c) and (d) agree 
well, e.g. at xu p /c = 1 and yu p /c = 0.5 or at xu p /c = 0.5 and yuj p /c = 4.5. We find 
that \Ej_\ w l-BjJ/lO in various locations and the electric force on an electron moving 
with the speed c/10 will equal the magnetic one. We decompose the total force into 
its constituents in Fig. 5(e,f). The major force contribution arises from the magnetic 
tension force. However, the omission of the magnetic pressure gradient in Fig. 5(e) 
alters the force distribution and results in clear differences with regard to \E±\. It is 
thus obvious that both components provide important contributions to the total electric 
field. No clear connection between B\\ and the other components is visible in Fig. 5. 

The currents can provide further information about the coupling between B\\ and 
the other field components. We subdivide for this purpose the current into a Jii = J z and 
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J± — Jx + iJy The power spectrum of both components is calculated by a 2D Fourier 
transform over the x—y plane followed by the azimuthal integration and a normalization 
by the maximum power of Jn . This provides us in analogy to the field distribution Fig. 
3 with the time-dependent power spectra of Jii and Jj_, which we display in Fig. 6. 




(a)log 10 kcte p 




Figure 6. (Colour online) The power spectra P(k, t) of the current components Jy 
(a) and J_l (b). Both spectra are normalized to the maximum power of Jn. The 
wavenumber, the time and the colour scale are 10-logarithmic. Ovcrplotted are two 
curves k oc t , which confine the wavenumber interval with the maximum power. The 
curves are the same in both panels. 

The expectation that Jy ^> J± due to \B±\ ^> \B«\ is confirmed. Both current 
distributions reveal that the wavenumber intervals, in which the power spectra of the 
currents peak, can be confined by two curves that follow k oc t~ x . The exact dependence 
of the wavenumbers, at which the power spectra of the currents peak is, however, not 
a power of t. The latter would correspond to a straight line in the double logarithmic 
diagram. Both current distributions rapidly expand to lower k after tu p ~ 125. This 
broadening of the current distribution is probably connected to that of the fields at late 
times in Fig. 3. A finite box effect may be responsible for this sudden broadening, 
because the largest observed wavelengths are not longer small compared to the box 
size. We discuss one such finite box effect below. However, this broadening sets in well 
after the saturation time tu p ~ 75 and these two processes are not connected. The 
signal-to-noise ratio of the currents is sufficient to reveal the cut-off of both spectra at 
logio(kc/u p ) = log 10 v^4, which is in line with Fig. 1. 

The overplotted curves k oc t^ 1 show that the wavenumber intervals, in which 
the power spectra of both current components peak, are practically identical. This 
confirms their connection. The abrupt growth of the power spectrum of J± at tu p ~ 70 
furthermore suggests that it is caused by a non-linear process. The immobile ions imply, 
that the Jj_ is caused by the bulk motion of the electrons in the simulation plane. The 
electrons can be accelerated either by E± or by their deflection by B± from the parallel 
direction into the perpendicular plane. The electromagnetic fields that accelerate the 
electrons in the perpendicular plane are tied to the spatially non-uniform Jn, which 
is driven by TAWI. The fields' extent is comparable to the spatial size of the current 
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filaments. It is thus not surprising that the current structures in J± have a typical size 
that is comparable to that of the structures in J». 

The time evolution of the current distributions J\\(x,y) (upper panel) and J±(x,y) 
(lower panel) are shown in the movie 2 for a subsection of the simulation box. The 
merging of filaments in Jy can be observed. Some merge to chains. If stable chains form, 
that spread across the entire simulation box with its periodic boundary conditions, then 
this would result in finite box effects. The formation of such long chains of filaments may 
cause the sudden spread of the currents to low k in Fig. 6 and in the electromagnetic 
fields in Fig. 5. The current structures in J± are driven by the fields and the current 
is dissipated away, resulting in their limited lifetime. The structures in J±, e.g. the 
current eddies in the lower panel of movie 2, give rise to the B\\ in Fig. 5. 

3.4- The velocity distribution 

The strong Jy should visibly modulate the phase space distribution f(x,y,v z ). It is, 
however, unclear if the phase space density remains compact or if beams form, i.e. 
distinct electron distributions for a given position along x and y. Figure 7 displays the 
phase space distribution as a function of v»(y) for a fixed xu p /c = 16.8. This phase 
space diagram is shown for the four simulation times tu p = 118, 198, 278 and 358. The 
velocity distribution, which has initially been spatially uniform, has been transformed 
into one, for which the mean speed varies as a function of the position. This velocity 
rearrangement has been accompanied by an energy transfer into the magnetic field. 

At tujp = 118, just after the TAWI has saturated, the velocity distribution oscillates 
in space with a high frequency, which corresponds to high k in Fig. 6. The spatial 
oscillation frequency decreases with time, while the thermal spread of the electron beam 
increases. The 2D distribution is scanned along the x-direction shown in Movie 3 at 
tu p = 358, showing that the velocity oscillations along the slice xu p /c = 16.8 in Fig. 7 
are representative for those in the simulation box. The actual velocity spread, within 
which we find significant numbers of electrons (note the linear density scale) is the same 
during the non-linear phase with i\v z \/v t h\\ % 1-5). Faster particles do exist. Initially 
the Maxwellian along the parallel direction could be populated with the given statistical 
plasma representation of 160 particles per cell up to \v max \ ~ 2.5i> t fcy. The low number 
density of the fast particles implies, however, that they do not contribute much to Jy. 
The simulation furthermore evidences, that no electrons are accelerated to high speeds 
during the non-linear evolution of the TAWI despite the large initial value of A. 

The increase with time of the wave length of the spatial oscillations in Fig. 7 
is related to a merging of the filaments in the Jy distribution (Movie 2), which could 
not be observed to this extent in the ID simulation of our previous work [17]. In 
the latter work, the particles have instead been pushed together by the magnetic 
pressure gradient, yielding a layered structure in the velocity distributions. In our 2D 
simulation the filament merging destroys this effect. In agreement with these results, the 
structure becomes more diffuse and the normalised current density j z n e /q cp c = J z /q cp c, 
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Figure 7. (Colour online) The colour plot shows the normalised velocity distribution 
v z{y)/ v th\\ f° r xcjp/c = 16.8 at the simulation times tuj p = 118 (a), 198 (b), 278 (c) 
and 358 (d). Overplotted is the normalised current density j z {y)n e /q cp c. 

overplotted in Fig. 7, with the current follows the structure of the velocity distribution. 
4. Discussion 

The thermal anisotropy-driven Weibel instability (TAWI), which is an important seed 
mechanism for magnetic fields in astrophysical environments and in laser-generated 
plasmas [1,2,10], has been investigated here with the help of a particle-in-cell (PIC) 
simulation for mobile electrons and for immobile ions. The plasma has initially 
been spatially uniform and we have set the initial electric and magnetic fields to 
zero. The thermal speed v t u\ of the electrons along one direction has been 20 times 
larger than the equivalent v t h± in the perpendicular plane. The anisotropy parameter 
A = v th\\/ v th± ~ 1 = 399 is thus very large. It allows us to examine the field growth and 
the electron thermalization under extreme conditions. It also provides us with a good 
signal-to-noise ratio of the electromagnetic fields in simulation. Thus, the results are 
valid for high anisotropies only. We plan to discuss the low anisotropy case in future 
work. 

The main purpose of our work has been to determine the source mechanism of the 
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electric field during the non-linear stage of the instability in more than one dimension. 
This electric field has recently received attention [15, 16], because the electric force is 
not small compared to the magnetic one. Both field components are thus important 
for the saturation of the TAWI, while typically only the magnetic field is considered 
[11, 23]. We have previously determined that the electric field is driven by the magnetic 
pressure gradient force, if the wave spectrum is limited to one dimension [17]. This force 
component is dominant for the filamentation instability in ID and 2D simulations [18], 
but it was not clear if this finding holds also for the TAWI. 

We summarize our results: We have compared the interval of unstable wave 
numbers in the simulation with the corresponding solution of the linear dispersion 
relation. Both agree prior to the saturation of the TAWI. Thereafter, the merging of the 
current filaments implies that the peak in the power spectrum of the current moves to 
lower wave numbers. This characteristic wave number decreases approximately linearly 
with an increasing time. The scale size of the filaments and the coherence length of 
the magnetic and the electric fields thus increases approximately linearly in time in the 
position space, until a sudden broadening sets in, which we have attributed to finite box 
effects. This evolution of the filament size contrasts the results of ID simulations, where 
mergers are possible only until the magnetic field becomes strong enough to keep the 
filaments with oppositely directed currents separated [11, 15, 17]. The repelling filaments 
can go around each other in a 2D plane and they continue to merge with other filaments, 
which have the same direction of the current vector. 

We have found that in a 2D simulation the magnetic tension force becomes stronger 
than the magnetic pressure gradient force, which clearly distinguishes this non-linear 
system from that driven by the beam filamentation instability. Both forces arise from the 
interaction of the net current driven by the TAWI with the magnetic field it generates. 
The electric field strength and distribution in the simulation plane resemble that of the 
force distribution obtained from the summation of the magnetic tension and pressure 
gradient force, at least in what concerns the strongest electric field structures. This 
electric field reaches an amplitude that makes it equally important for the dynamics 
of the slow electrons as the magnetic field. Its energy density remains, however, well 
below the magnetic one, which is in agreement with previous simulations. We find 
that the estimate for the magnetic field amplitude, which lets the TAWI saturate, that 
is based on the magnetic trapping mechanism is still a reasonable approximation for 
that observed in our simulation. This has also been reported for ID simulations of the 
non-relativistic TAWI [15] and for the filamentation instability [18]. 

The electromagnetic fields generated by the TAWI drive electron currents in the 
simulation plane. These currents result in the growth of a magnetic field component, 
which is orthogonal to the simulation plane and parallel to the direction, along which 
the electrons are hottest. The growth of such a field component is not predicted by the 
solution of the linear dispersion relation and it is thus a purely non-linear process. This 
magnetic field component does also not grow in ID PIC simulations [17]. Its amplitude 
is lower than that in the simulation plane but it is large enough to result in complicated 
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3D magnetic patterns in a 3D simulation [20]. 

The ratio of the energy density of the two magnetic field components in our 
simulation plane and the total electron thermal energy density exceeds the expected 
limiting value of 1/12 [12]. It remains, however, just below twice that value. Such a 
peak energy density of the magnetic field in the simulation plane is reasonable, if the 
limit 1/12 applies to each magnetic degree of freedom. One magnetic component is 
considered in Ref. [12], while two components grow in our simulation plane. 
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